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EXECUTIVE  SUMMARV 


This  paper  describes  the  methodology  and  results  of  a 
project  to  determine  the  magnitude  and  amplitude  dependence  of 
the  errors  inherent  in  the  full-Doppler  mode  NAVSPASUR  phase 
difference  data.  The  purpose  of  this  work  was  to  generate  an 
accurate  model  of  these  errors  for  subsequent  use  in  simulation 
programs  designed  to  predict  the  performance  capabilities  of  both 
the  current  NAVSPASUR  system  and  various  proposed  upgrades.  We 
report  here  on  the  preliminary  results  obtained  from  an  analysis 
of  full-Doppler  mode  data  collected  by  the  NAVSPASUR  system 
during  one  four  hour  period.  We  plan  in  the  near  future  to 
process  a  much  more  comprehensive  data  base  in  order  to  refine 
the  error  model  and  to  allow  us  to  address  the  possible  sources, 
of  these  errors.  No  attempt  has  been  made  in  the  current  work  to 
identify  the  source  of  the  errors  or  to  address  whether  (and/or 
which)  operational  improvements  could  be  made  to  reduce  their 
magnitude. 

The  main  results  of  this  project,  described  in  detail  in 
this  report,  are  as  follows: 

1.  For  a  given  satellite  pass  data  set  (hereafter  scan) 
the  phase  difference  errors  for  each  antenna  array  have 
both  random  and  systematic  components.  For  signals 
which  are  more  than  about  6  dB  above  the  noise  level 
the  systematic  errors  are  about  50%  larger  than  the 
purely  random  errors.  These  systematic  phase  difference 
offsets  are  relatively  constant  from  scan  to  scan,  at 
least  for  time  intervals  up  to  4  hours.  Their  behavior 
mimics  (and  they  may  in  fact  be)  phase  calibration 
errors.  For  strong  signals  the  random  errors  are 
comparable  in  size  to  the  phase  data  quantization  size 
(5.6  degrees).  However,  since  systematic  errors 
dominate  for  high  signal  to  noise,  the  phase 
quantization  does  not  appear  to  limit  overall  accuracy 
at  this  time. 

2.  The  random  phase  difference  errors  appear  to  be 
independent  of  actual  or  projected  baseline  length  or 
orientation  (N-S  vs.  E-W) .  The  systematic  errors  show 
a  bias  towards  larger  errors  on  the  shorter  baselines. 

The  net  (both  systematic  and  random)  phase  difference 
errors  are  normally  distributed  with  a  dispersion 
characterized  by  a  single  parameter  -  the  received 
signal  strength.  No  significant  station-to-station 
variations  are  apparent,  and  there  do  not  appear  to  be 
any  temporal  (day/night)  variations  of  significance. 
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The  RMS  phase  difference  error  as  a  function  of 
received  signal  strength  is  well-behaved,  and  its 
dependence  has  been  modelled  with  a  fourth-order 
polynomial  (see  Figures  2  and  3) . 

3.  An  effective  noise  floor  of  approximately  -152  dBm  is 
confirmed  as  representative  of  the  current  system. 
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I.  System  Description  and  Analytic  Approach 


The  fundamental  NAVSPASUR  measured  quantity  is  the  received 
signal  phase  at  each  antenna  array  at  a  given  receiver  site.  The 
raw  phases  from  each  of  the  n  antenna  arrays  at  a  receiver  site 
are  combined  to  produce  n-1  phase  differences,  each  taken  with 
respect  to  the  designated  reference  antenna.  The  deployment  of 
antenna  arrays  along  each  of  two  orthogonal  directions,  together 
with  the  use  of  antenna  pairs  (hereafter  baselines)  with 
different  spatial  separations  (baseline  lengths) ,  allows  for  the 
simultaneous  determination  of  two  direction  cosines  (cos  a,  where 
Ok  is  the  angle  between  the  baseline  vector  and  the  satellite 
position  vector  in  station  coordinates)  at  each  receiving  site. 
Further,  time-sequence  observations  of  the  phase  differences 
u^bociated  with  each  baseline  allow  for  the  calculation  of  two 
direction  cosine  rates,  d(cos  a)/dt,  which  are  functionally 
related  to  the  satellite's  velocity. 

Because  the  direction  angles  and  angle  rates  are  derived 
directly  from  the  phase  difference  data  rather  than  the  individual 
antenna  phases,  we  focus  here  on  determining  the  errors  on  the 
phase  differences  as  opposed  to  the  raw  phases.  The  RMS  phase 
difference  errors  due  to  random  raw  phase  errors  are  expected  to 
be  42  larger  than  the  RMS  raw  phase  errors.  Further,  the  phase 
differences  will  also  include  systematic  errors  due  to  phase 
biases,  or  offsets,  which  would  not  appear  in  an  analysis  of  raw 
phase  data. 

The  phase  difference  data  from  a  given  receiver  site  allows 
for  the  calculation  of  four  parameters  which  provide  information 
on  the  satellite  position  and  velocity.  In  order  to  fully 
specify  the  position  and  velocity  of  the  satellite  at  some  given 
time  tQ/  a  total  of  six  parameters  are  required.  The  NAVSPASUR 
system  meets  this  requirement  by  combining  measurements  of  four 
parameters  from  each  of  two  or  more  spatially  separated  receiving 
sites.  It  is  important  to  note  here,  however,  that  because  of 
the  particular  geometry  of  the  NAVSPASUR  system,  the  parameters 
as  determined  from  the  different  receiving  stations  are  to 
varying  degrees  linearly  dependent  on  one  another,  so  that  the 
satellite  position  and  velocity  are  not  necessarily  well 
determined  in  all  cases.  For  this  reason  the  relationship 
between  single  station  performance  and  system  performance  is  a 
complex  issue,  and  its  treatment  is  deferred  to  a  subsequent 
report,  which  will  build  on  the  results  contained  herein. 

The  phase  difference  data  from  which  the  direction  cosines 
and  cosine  rates  are  obtained  are  contaminated  by  random  noise 
and  systematic  measurement  and  calibration  errors.  In  the 
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absence  of  any  noise  or  systematic  errors,  the  direction  cosines 
and  cosine  rates  are  overdetermined  by  the  available  data.1  In 
the  presence  of  noise  and  systematic  errors,  the  problem  becomes 
one  of  obtaining  a  best  estimate  of  the  desired  parameters  given 
the  available  data.  While  a  variety  of  estimation  techniques 
exist,  perhaps  the  most  common  and  straightforward  technique  is 
the  method  of  least  squares.  With  this  method  one  adopts  as  the 
best  estimators  those  values  of  the  independent  parameters  which 
minimize  the  sum  of  the  squares  of  the  residual  errors.  A 
substantial  body  of  mathematics  exists  in  support  of  the  least 
squares  method,  and  it  is  extensively  discussed  in  many 
statistical  textbooks.  For  completeness,  a  description  of  the 
methodology  is  included  here,  following  the  discussions  contained 
m  Lewis  and  Odell  (1971)  and  Mathews  and  Walker  (1970). 

II.  The  Method  of  Least  Squares 


a).  Basic  formulation  of  the  nonlinear  least  squares  problem 

The  method  of  least  squares  provides  a  means  of  obtaining  a 
best  estimate  of  a  set  of  parameters,  represented  by  a  parameter 
vector  *  of  dimension  m,  given  an  associated  data  vector  y  of 
dimension  n  >  m  such  that  the  vector  /(*)  gives  an  approximate 
correspondence  between  *  and  y.  In  general,  the  vector  /  is  such 
that  the  parameters  *  may  be  expressed  nonlinearly  in  /.  Where 
the  functional  form  of  /  is  known,  the  problem  reduces  to  one  of 
determining  the  parameter  vector  i  which  most  closely  describes 
the  n  relationships  expressed  by  the  vector  equation  y  =  f(i) ,  in 
the  sense  that  the  sura  of  the  squares  of  the  residuals, 

Q (*)  =  [y  -  /(*)]T[y  -  '(*)]  =  I  [ y ,  -  /,(*)]2  d) 

J=1  1  1 

is  minimized.  This  requirement  is  satisfied  when  the  partial 
derivatives  of  Q  with  respect  to  the  parameters  *,•  equal  zero,  i.e., 

Gt(*)  s  -  -2  [y  -  /(*)  ]T[*t-(*)  ]  =0,  i  =  1,  2,  .  .  .  ,  m.  (2) 

1.  For  each  receiving  station  there  are  *  determinations 
of  each  direction  cosine  and  *-l  determinations  of  each 
direction  cosine  rate  where,  typically,  20<  *  <55. 
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(3) 


Here,  $.(*)  is  an  n-dimens ional  vector  given  by 


A  _  ^(*) 
*i(*>  ~  a*; 


In  the  case  where  the  output  vector  y  is  expressed  linearly 
in  *,  the  conditions  given  by  equations  (2)  yield  m  linear 
equations  in  m  unknowns,  which  can  be  solved  directly  for  the  *■. 
When  the  function  /  is  nonlinear  in  *,  the  usual  approach  to 
minimizing  equation  (1)  is  to  linearize  equations  (2)  by  choosing 
a  suitable  first  estimate  *0  of  *  and  replacing  /(*)  in  equation 
(1)  with  its  first-order  Taylor  series  expansion  about  *0, 


/(*)  =  /»(*)  =  /(*0)  *  S  (**  -  **>♦*(♦„) 

k- I 


Substituting  /?(*)  for  /(*)  in  equation  (1)  and  setting  the 
partial  derivatives  0t(*)  equal  to  zero  then  yields  the  following 
system  of  linear  equations: 


*1*2 


<t>T  $  $ 

*irr  l  *nr2 


<bT<b 
"  m 


A 


where  <t>i  =  and  Qi  =  Gj(*0)  .  Equation  (5)  can  be  written  as 

the  matrix  equation 

F  d  =  4  Q  (6) 

with  the  solution 


of  =  -i  F  'O 
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The  elements  dt  =  -  *0l  )  ,  i  =  1,  2,  .  .  .  ,  m  represent 
corrections  to  the  first  estimate  *0  of  *.  Equation  (7)  may  be 
solved  for  the  vector  d  using  standard  matrix  inversion 
techniques.  From  the  corrections  d,  one  obtains  an  improved 
estimate  *,=*0  +  0?  for  the  parameter  vector  *.  Replacing  *0  in 
equation  (4)  with  the  above  process  is  repeated  to  obtain  a 
new  estimate  i2  f or  *•  This  iterative  procedure  continues  until 
one  obtains  sufficiently  accurate  estimates  for  each  of  the 
parameters  *t-. 


b).  Weighting  data  of  nonuniform  precision 

In  the  above  analysis  we  have  assumed  that  all  data  points 
have  the  same  relative  precision,  so  that  all  points  contribute 
equally  to  obtaining  the  least  squares  solution.  When  the  data 
are  not  all  of  equal  precision  the  least  squares  best  estimate  is 
obtained  when  the  weighted  sum  of  the  squares  of  residuals,  i.e. , 

<?(*)  =  -  ',•(*>  ]2  (8) 
i=l  1  1  1 

are  minimized.  Here  the  w-  represent  the  weights  for  the 

individual  data  points.  If  the  errors  for  all  measurements  of  a 
given  precision  are  assumed  to  be  normally  distributed,  then  the 
relative  weight  of  those  data  points  is  inversely  proportional  to 
their  variance,  so  that  Wj  =  l/o;2.  If  the  variances  0?  are  not 

known  a  priori,  one  may  proceed  by  making  a  reasonable  first 
estimate  of  their  values,  perform  the  least  squares  fitting, 
calculate  refined  values  of  the  variances  from  the  results,  and 
use  these  newly  obtained  values  as  input  to  a  second  running  of 
the  least  squares  routine. 

c).  Error  estimation  and  correlation  of  errors 

Once  a  least  squares  estimate  of  the  parameter  vector  *  has 
been  obtained,  it  is  useful  to  investigate  the  probable  errors  in 
the  determined  parameters.  If  the  experiment  under  investigation 
were  repeated  a  large  number  of  times  in  exactly  the  same  fashion 
and  with  the  same  error  distribution,  we  could  calculate  the 
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mean-square  deviation  of  i  for  any  individual  experiment  from  the 
grand  average  ?  taken  over  all  experiments.  To  find  the  errors 
in  *,  then,  we  need  to  calculate  the  mean-square  deviations 
<(*t  -  Since  (*t  -  *•)  =  ( di  -  d()  ,  we  have,  from  equations 

(7)  ,  (8)  ,  and  (3)  , 

m  -  m  fi 

<*(  -  *<)  -  S  <r‘>tt  (o*  -  o*)  -  11  (/"'I*  %  (y,  -  v,)  .  (»> 

where  ^  Since  the  individual  measurements  y^.  are 

assumed  to  be  statistically  independent,  the  expectation  values 
of  the  errors  must  satisfy 


<  - 


(yn  -  y. 


r>)  = 


6*  = 


(10) 


so  that,  in  general, 


_  _  y  772  /l  /7l  /l  _  *1 

*<  -  *<><*/  -  W  =  222,2  ('-'>*  <'%  ♦*,  V  -y  -7  ij  */, 


“  22  <r'>*  V 


=  (r')„  . 


The  standard  error  A*,  in  is  therefore 


(ID 


A* 


1/2  - - - 

i  '  <(*i  -  5i>2)  -  Nl(r')„  • 


(12) 
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The  matrix  F'{  represents  an  error  matrix  whose  diagonal 
elements  are  the  squares  of  the  standard  errors  in  the  parameters 
*£.  The  fact  that  the  cross  terms  -  *,)>,  with  i  f  l 

are  not  in  general  equal  to  zero  implies  that  the  parameters 
are  not  statistically  independent  and  that  the  errors  are 
correlated.  The  expected  error  in  any  quantity  which  is  a 
function  of  the  parameter  vector  i  may  be  determined  from  the 
error  matrix.  For  a  given  function  g (i) ,  the  standard  error  in 
g(i)  due  to  the  uncertainties  in  the  parameters  *  may  be  found 

from  the  relation  <[£(*)  -  g( *)  ]  /  =  <[£(*  -  *)  ]  /  by  direct 

substitution. 

III.  The  Application  of  Least  Squares  to  the  NAVSPASUR  System 

We  have  applied  a  nonlinear  least  squares  analysis  with 
nonuniform  data  weighting  to  full-Doppler  mode  NAVSPASUR  phase4 
data  in  order  to  determine  the  distribution  of  the  residual  phase 
difference  errors.  For  each  data  scan  (defined  as  the  data 
collected  by  one  receiving  station  from  a  single  pass  of  a  target 
satellite  illuminated  by  a  single  transmitter)  we  have  obtained  a 
best  fit  to  the  E-W  and  N-S  direction  angles  (hereafter  6  and  #, 
respectively)  and  direction  angle  rates  (hereafter  6  and  $)  at 
the  epoch  of  fence  crossing.  From  these,  we  calculate  the 
expected  value  of  the  differential  phase  for  each  datum  and 
derive  the  residual  phase  difference  errors.  We  then  analyze 
these  residual  errors  to  obtain  their  dependence  on  received 
signal  amplitude.  The  details  of  this  process  are  outlined  in 
the  following  sections. 

a).  Coordinate  system  and  basic  equations 

For  the  single-station  fits  to  the  data,  the  most  natural 
coordinate  system  in  which  to  work  is  the  local  station 
coordinates,  where  the  reference  antenna  is  at  the  origin  and  the 
orthogonal  axes  are  local  "north"  and  "east",  both  tangent  to  the 
geoid,  and  local  height,  or  up,  which  is  perpendicular  to  the 
surface  of  the  geoid.  Here  "east"  is  eastward  in  the  plane  of 
the  NAVSPASUR  great  circle,  and  "north"  is  normal  to  the  great 
circle  in  a  northward  direction.  Figure  1  depicts  the 
orientation  of  a  local  station  coordinate  system  on  the  earth's 
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surface.  For  this  coordinate  system,  all  antenna  offsets  are  to 
the  "east"  and/or  "north",  with  no  vertical  component  (i.e. ,  the 
antennas  all  lie  in  the  tangent  plane)  .  The  offsets  from  the 
coordinate  origin  can  be  obtained  directly  from  the  station 
survey  plats.  Further,  the  station  location  in  geocentric 
coordinates  and  the  rotation  vectors  needed  to  transform  the 
satellite  geocentric  position  and  velocity  to  local  station 
coordinates  are  known  and  can  be  checked  for  consistency  against 
those  currently  used  by  NAVSPASUR.  In  this  coordinate  system  the 
direction  angles  8  and  #  and  the  distance  to  the  satellite  rsat  can 
be  expressed  as 


9 ( t )  =  tan" 


rE(t) 

rH(t) 


*(t)  =  sin 


rN(t) 

rsa«(t) 


(13) 


rsat(t)  =  (  rN(t)2  +  rE(t)2  +  rH(t)2  )w2 


where 


(t-t0)2  (N  •  R) 

(t-t0)2  (E  •  R)  (14) 

(t-t0) 2  (H  •  R)  . 

Here  t0  is  the  epoch  of  fence  crossing,  rz  and  v2  are  the 
satellite  position  and  velocity  in  local  station  coordinates,  G 
is  the  gravitational  constant,  Mg  is  the  mass  of  the  earth,  R  is 
the  satellite's  position  in  geocentric  coordinates,  and  the  N,  E, 
H,  are  the  unit  direction  vectors  of  the  three  local  station 
coordinate  axes.  The  angles  used  here,  together  with  the 


rN(t)  =  rN(tD)  +  vN(t0)(t-t0)  +  ^ 


rE(t)  =  rE(t0)  +  vE(t0)  ( t-t0)  +  ^ 


rH(fc)  ~  r„(t0)  +  vH(t0)(t-t0)  + 


7 


satellite  range,  represent  an  orthogonal  polar  coordinate  system 
ir  which  9  measures  the  rotation  angle  in  the  great  circle  plane 
and  #  measures  the  angle  normal  to  the  great  circle.  As  such, 
the  angles  9  and  #  do  not  exactly  correspond  to  the  angles 
measured  by  the  two  direction  cosines  at  the  station  location. 
However,  for  the  particular  geometry  of  the  NAVSPASUR  system, 
where  the  satellite  is  illuminated  only  within  a  few  arcminutes 
of  the  great  circle  plane,  these  angles  do  closely  correspond  to 
those  measured  by  the  direction  cosines.  For  computational 
simplicity  it  is  preferable  to  work  with  these  angles  as  opposed 
to  working  directly  with  the  direction  cosine  angles.  This 
choice  of  angle  definitions  has  no  effect  on  the  results  of  our 
analysis . 

Over  the  short  duration  of  the  satellite's  illumination  by 
the  transmitter  beam  (<  1  sec),  the  last  (acceleration)  terms  in 
equations  (14)  can  safely  be  ignored.  The  direction  angle  rates 
and  range  rate  may  be  obtained  by  taking  the  time  derivatives  of 
equations  (13),  where  again  the  acceleration  terms  can  safely  be 
ignored,  so  that  we  have 

A/_  VE  (t0)  cos{9  (t )  )  -  vH(t0)sin(9(t)  ) 

0(t) - Rtj 


*(t)  = 


vN(tc)-vE(t0)sin{9(t)  )  sin{  $  ( t)  )  -vH  ( t0)  cos  { 9  ( t)  )sin(*(t)  ) 

P(t) 


(15) 


*s.t(t> 


vN(t0)rN(t)  +  vE(t0)rE(t)  +  vH  ( t0)  rH  (t) 

rsat  ( t) 


where 


p(t)  =  [  rH  (t)  +  rE  (t)  ]1/2  . 


For  the  least  squares  fitting  procedure,  equations  (13)  and 
(15)  provide  a  means  of  obtaining  a  first  estimate  of  the  angles 
and  angle  rates  from  the  predicted  satellite  position  and 
velocity  at  the  epoch  of  fence  crossing.  Because  the  single 
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station  data  can  effectively  constrain  only  these  four 
parameters,  we  must  fix  the  two  remaining  parameters,  rsa,  and  r  . 
The  logical  choice  for  these  parameters  is  to  fix  them  to  their 
a  priori  values  based  on  the  ephemeris  prediction.  Since  the 
phase  differences  are  insensitive  to  range  and  range  rate,  this 
approximation  introduces  no  significant  error  into  our 
calculations.  With  this  assumption,  then,  the  satellite  position 
at  time  t  on  subsequent  iterations  is  given  in  terms  of  the 
current  estimates  of  the  direction  angles  and  direction  angle 
rates  by 


Mt)  =  sin{*n(t)  } 

rE  ( t)  =  rsat(t)  cos{*n(t) }  sin{6n(t)  )  (16) 

rH(fc)  =  rsai  ( 1 )  cos{*n(t)  )  cos(0n(t)  ) 

where 

=  rsa.(V  +  *s„(t0)(t  ~  fco) 

*«(t)  =  *«(t0)  +  *„(t0)  (t  -  t0)  (17) 

e„(t)  =  en(t0)  +  e„(t0)  (t  -  t0)  . 

Here,  the  subscript  n  has  been  added  to  the  angles  and  angle 
rates  to  indicate  their  value  after  the  nth  iteration.  Equations 
(16)  and  (17)  are  the  basic  equations  from  which  the  satellite 
position  (and  the  expected  phase  differences)  are  calculated  for 
all  least  squares  iterations  after  the  first. 

b).  Derivation  of  the  ’’ideal”  phase  difference 

For  a  satellite  of  known  position  and  velocity,  the  phase 
difference  (in  rotations)  between  two  receiving  antennas  at  known 
positions  is  simply  the  difference  in  the  respective  path  lengths 
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from  the  satellite,  expressed  in  wavelengths.  For  the  satellite- 
receiver  distances  typical  for  the  FAVSPASUR  system,  one  cannot 
use  a  plane  wave  approximation  over  the  longest  baselines,  and 
the  respective  path  lengths  must  be  explicitly  calculated.  Since 
the  reference  antenna  has  been  taken  as  the  origin  of  the  local 
station  coordinate  system  used,  this  phase  difference  is  simply 


■  sat 


-  rC  -  Ir 


sat1 


(18) 


where  is  the  phase  difference  at  the  fcth  antenna  relative  to 

the  reference  antenna,  rk  is  the  position  vector  of  the  *th 
antenna,  and  krec  is  the  wavelength  of  the  received  ( i.e. ,  Doppler 
shifted)  signal.  In  actuality,  the  NAVSPASUR  system  only 
measures  antenna  phases  modulo  1  (i.e.,  0  to  1  rotations,  or  o  to 
360  degrees) .  In  the  least  squares  fitting  program  we  perform 
the  appropriate  modular  arithmetic  and  express  all  phase 
differences,  both  observed  and  ideal,  in  the  range  -0.5  to  +0.5 
rotations . 

At  the  beginning  of  each  iteration  of  the  least  squares 
routine  we  calculate  the  ideal,  or  expected,  phase  differences 
based  on  the  satellite  position  as  derived  from  equations  (16) 
and  (17) . 


c).  Calculation  of  the  phase  difference  partial  derivatives 


In  addition  to  the  the  expected  phase  differences,  we  also 
need  for  the  least  squares  routine  the  partial  derivatives  of  the 
phase  differences  with  respect  to  the  parameters  to  be  solved  for 

(#,  0,  4,  and  0).  We  obtain  these  using  a  finite  difference 
scheme,  where,  for  example,  the  partial  derivative  with  respect 
to  *  of  the  phase  difference  for  the  baseline  associated  with 
the  *th  antenna  is  given  by 


where  we  take  6 i  =  i/10s.  The  partials  for  the  three  remaining 
parameters  are  calculated  in  an  analagous  manner. 

d) .  Data  calibration 

For  each  data  tape  received  from  Dahlgren  we  extract  the 
half-hourly  phase  calibration  data  scans.  For  each  calibration 
scan  we  compute  the  differenced  phases  with  respect  to  the 
reference  antenna  for  each  antenna  and  timeline,  and  average  all 
timelines.  Prior  to  the  averaging  process,  the  data  are  screened 
to  remove  any  spurious  data  points  which  occasionally  occur  at 
the  end  of  the  cal  scans.  Given  the  good  phase  stability  of  the 
system,  these  points  are  easily  identified  and  removed.  If  the 
data  tape  spans  more  than  one  calibration  scan  for  a  given 
receiver  site,  the  differenced  phases  for  each  baseline  are 
averaged  across  all  available  calibration  scans,  which  again  is 
possible  owing  to  the  excellent  phase  stability  of  the  system. 

Once  the  phase  calibration  values  for  all  the  receiving 
stations  have  been  established,  these  calibration  values  are 
subtracted  from  the  individual  satellite  data  scans.  Finally, 
the  phase  differences  for  each  baseline  are  computed  for  each 
data  scan.  Only  data  scans  for  which  the  target  satellite  has 
been  identified  as  a  known  catalog  object  are  retained  for 
processing,  since  the  nonlinear  least  squares  analysis  requires 
reasonably  good  initial  estimates  of  the  position  and  velocity  to 
achieve  convergence.  The  calibrated  phase  difference  data, 
together  with  the  corresponding  received  signal  amplitude,  time 
of  observation,  Doppler  frequency,  catalog  position  and  velocity, 
and  epoch  of  fence  crossing  are  written  to  a  calibrated  data  file 
which  is  used  as  input  to  the  least  squares  fitting  program. 

e) .  Data  weighting 

Since  the  purpose  of  this  project  was  to  derive  an  accurate 
error  model  for  the  full -Doppler  mode  NAVSPASUR  data,  we  did  not 
have  an  initial  error  model  available  to  use  in  the  data 
weighting.  From  an  examination  of  the  calibrated  data  scans  we 
produced  an  initial  rough  estimate  of  the  error  model,  fitting 
the  estimated  errors  to  a  piecewise  linear  function  of  signal 
strength.  This  initial  model  was  used  as  the  error  weighting 
model  for  the  first  pass  least  squares  fitting  of  the  data,  and 
from  the  first  pass  we  generated  a  fourth-order  polynomial  error 
model.  This  new  model  was  then  used  as  the  input  error  model  for 
data  weighting  in  a  subsequent  pass  of  the  least  squares  routine. 
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In  order  to  test  the  sensitivity  of  the  results  to  the 
assumed  error  model  we  also  ran  the  data  base  with  no  (i.e., 
uniform)  weighting  of  the  input  data.  A  comparision  of  the  final 
error  models  for  each  of  these  three  cases  demonstrated  that  the 
form  of  the  final  error  model  was  not  very  sensitive  to  the  form 
of  the  input  error  model  used  for  the  data  weighting. 

/).  Least  squares  fitting  to  the  data 

Initial  estimates  of  the  angles,  angle  rates,  range,  and 
range  rate  are  calculated  from  the  a  priori  (catalog)  satellite 
position  and  velocity  contained  in  the  data  scan  header  using  the 
relationships  given  in  equations  (13)  and  (15) .  In  the  least 
squares  fitting  process  we  fix  the  range  and  range  rate  to  their 
initial  values,  since  data  from  a  single  receiving  station  can 
constrain  only  the  remaining  four  parameters.  The  angles  and 
angle  rates  form  a  four  dimensional  parameter  vector  *  which  we 
wish  to  estimate  in  the  least  squares  sense.  From  the  initial 
estimate  *Q  of  *  we  calculate  the  expected  phase  differences  p(*0) 
and  the  partial  derivatives  of  the  phase  differences  $(*)  .  For 
each  datum  we  also  calculate  the  data  weight  w  =  1/a2  based  on  a 

t 

piecewise  linear  initial  estimate  of  the  error  model.  We  then 
construct  the  normal  equations  as  described  in  SII,  yielding  four 
equations  in  the  four  unknown  parameters,  i.e. ,  the  four 
corrections  dt  to  the  initial  estimates  *0.  From  these  we  compute 
new  estimates  which  are  used  as  the  starting  point  for  the 
second  iteration.  These  new  estimates  are  used  to  obtain  the 
satellite  position  via  equations  (16)  and  (17),  and  the  remainder 
of  the  process  continues  as  for  the  initial  iteration. 

We  continue  iterating  in  this  fashion  until  the  solution 
converges  or  until  a  specified  maximum  number  of  iterations  are 
completed.  The  convergence  criterion  we  adopted  was  that  the 
magnitude  of  each  of  the  corrections  dt  at  the  end  of  an  iteration 
be  less  than  or  equal  to  the  estimated  errors  in  the  corrections 
Adt  =  Atj-.  If  this  convergence  test  is  not  met  within  three 
iterations  the  fitting  process  is  abandoned  and  we  proceed  to  the 
next  data  scan.  The  three  iteration  limit  was  selected 
empirically  based  on  our  experience  in  generating  fits  to 
NAVSPASUR  data.  When  the  satellite  signal  is  sufficiently  above 
the  noise  level  and  the  identification  in  the  scan  header  is 
correct,  the  least  squares  routine  typically  converges  in  one  or 
two  iterations.  When  the  routine  is  unable  to  converge  within 
three  iterations,  our  experience  shows  that  it  is  unlikely  to 


converge  at  all.  The  failure  to  achieve  convergence  is  due  in 
almost  all  cases  to  poor  signal-to-noise .  Therefore,  to  prevent 
the  waste  of  excessive  amounts  of  computer  time,  we  set  the 
convergence  limit  to  three  iterations.  All  data  scans  for  which 
the  program  failed  to  reach  convergence  were  retained  for  future 
analysis.  The  non-convergent  scans  represent  approximately  8%  of 
the  total  number  of  scans  processed. 

Once  the  least  squares  fits  to  the  observational  data  were 
complete,  the  residual  phase  difference  errors  for  each  datum 
were  computed.  These  residuals,  together  with  the  corresponding 
received  signal  amplitudes  and  the  data  scan  headers,  were  saved 
to  disk  files  for  further  analysis. 

IV.  Analysis  of  the  Phase  Difference  Errors 

The  analysis  we  present  here  is  based  on  a  four  hour 
NAVSPASUR  phase4  data  tape  containing  all  observations  for  the 
period  20AUG87  01:50  UT  -  06:00  UT.  Since  the  nonlinear  least 
squares  fitting  routine  requires  a  good  initial  estimate  of  the 
satellite  position  and  velocity,  only  those  scans  for  which  a 
catalog  identification  was  assigned  were  processed.  This 
amounted  to  a  total  of  4,212  data  scans.  Prior  to  least-squares 
processing  of  the  data,  the  calibration  scans  from  the  data  tape 
were  carefully  checked  to  ensure  that  no  calibration  errors  were 
present  which  would  propagate  into  the  data  analysis.  The 
calibration  data  were  found  to  be  well-behaved  over  the  full  time 
interval  spanned  by  the  data. 

During  the  time  period  spanned  by  the  data  tape  the 
Hawkinsville  station  was  operating  in  the  '  J1 '  mode  instead  of 
its  normal  configuration.  Therefore,  we  have  excluded  the 
Hawkinsville  data  (925  scans)  from  our  analysis.  This  exclusion 
is  not  expected  to  significantly  affect  our  results,  since  we 
have  analyzed  Hawkinsville  data  from  other  shorter  time  period 
data  tapes  and  find  that  the  error  distribution  is  similar  to 
that  which  we  report  here  for  the  remaining  stations. 

Approximately  17%  of  the  data  scans  which  nominally  met  the 
convergence  criteria  produced  clearly  anomalous  results,  in  that 
the  mean  residuals  for  the  majority  of  the  baselines  were 
excessively  large  (typically  70  -  150  degrees),  while  the 
dispersion  in  the  residuals  was  normal  (typically  9  degrees) .  We 
believe  these  scans  to  be  cases  for  which  the  assigned  catalog 
identification  is  incorrect.  In  the  event  of  an  incorrect 
satellite  identification,  it  is  possible  (and,  in  fact,  likely) 
that  the  nonlinear  least  squares  routine  could  converge  on 
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incorrect  values  of  the  parameters.  This  occurs  because  the 
nonlinear  least  squares  method  requires  a  reasonably  good  first 
estimate  of  the  parameters  in  order  to  converge  on  the  correct 
answer.  Given  a  poor  first  estimate  of  the  parameters,  the 
algorithm  can  converge  on  a  local,  rather  than  global  minimum  of 
the  sum  of  the  squares  of  the  residuals. 

We  have  examined  several  of  the  anomalous  scans  and  have 
been  able  to  verify  in  many  cases  that  the  originally  assigned 
identification  was  incorrect.  While  we  have  not  been  able  to 
verify  that  this  was  the  case  for  all  the  anomalous  scans,  it 
does  appear  that  they  may  all  be  due  to  the  least  squares  routine 
converging  on  the  wrong  minimum.  To  avoid  contamination  of  our 
results  be  these  anomalous  scans,  we  have  eliminated  them  from 
our  data  analysis. 

After  filtering  the  data  to  remove  the  nonconvergent  and 
anomalous  scans,  and  all  scans  for  Hawkinsville ,  our  data  set 
consisted  of  3,055  scans.  We  have  analyzed  the  remaining  data  at 
three  levels  and  intercompared  our  results  from  each  level,  as 
described  below. 

First,  we  examine  the  data  on  an  individual  scan  basis.  At 
this  level  of  examination,  there  is  insufficient  data  to 
determine  the  amplitude  dependence  of  the  errors.  We  therefore 
have  combined  all  data  whose  amplitudes  exceed  the  noise  level  by 
more  than  6  dB  on  a  baseline-by-baseline  basis  and  perform  a 
statistical  analysis.  We  calculate  for  these  data  the  mean,  the 
dispersion  about  the  mean,  and  the  RMS  magnitude  of  the 
residuals.  If  the  noise  were  purely  random,  the  mean  of  the 
residuals  for  each  baseline  should  be  zero  within  the  statistical 
errors,  and  the  calculated  dispersions  and  RMS  deviations  should 
be  approximately  equal.  We  find  that,  in  general,  this  is  not 
the  case. 

Typically,  the  absolute  value  of  the  mean  of  the  residuals 
for  any  baseline  is  about  14  degrees,  and  the  dispersion  about 
the  mean  is  about  9  degrees.  The  deviation  of  the  means  from 
zero  is  statistically  significant  given  the  number  of  data 
points.  Further,  the  means  for  any  given  baseline  correlate  well 
from  scan  to  scan,  even  over  time  intervals  of  up  to  four  hours 
(the  longest  time  period  sampled  to  date) .  This  suggests  that 
there  are  in  fact  significant  systematic  effects  contributing  to 
the  phase  difference  errors,  a  point  we  will  return  to  below. 
The  magnitude  of  the  combined  systematic  and  random  errors  is 
significantly  larger  than  the  phase  quantization  size  in  the  raw 
data  (5?6) ,  suggesting  that  for  the  current  operation  phase 
quantization  effects  are  not  a  limiting  factor  to  the  system 
performance. 
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The  second  level  of  data  examination  occurs  at  the  station 
level,  where  data  for  the  whole  time  period  for  a  given  receiving 
station  is  combined  on  both  a  baseline-by-baseline  basis  and  on 
an  overall  station  basis.  At  this  level  we  are  able  to  bin  the 
data  by  received  signal  amplitude  and  determine  the  amplitude 
dependence  of  the  errors.  For  data  whose  amplitude  exceeds  the 
system  noise  level  by  more  than  6  dB,  we  find  that,  for  any  given 
baseline,  the  RMS  errors  are  constant  with  signal  strength  and 
are  typically  11-29  degrees.  The  significant  variations  we  see 
between  different  baselines  is  due  primarily  to  the  fact  that  the 
means  are  not  zero  within  statistical  errors,  presumably  due  to 
systematic  errors.  The  standard  deviation  about  the  mean, 
typically  9  degrees,  remains  fairly  constant  from  baseline  to 
baseline,  however. 

We  have  not  searched  in  detail  for  correlations  in  the 
apparent  systematic  errors,  but  a  visual  inspection  of  the  data 
indicates  that  the  magnitude  in  the  systematic  component  anti¬ 
correlates  with  baseline  length.  That  is,  the  shortest  baselines 
exhibit  the  largest  mean  deviation  from  zero.  The  interpretation 
of  this  trend  may  be  misleading,  however,  since  systematic  errors 
on  the  longest  baselines  can  be  compensated  for  by  relatively 
small  changes  to  the  best  estimates  of  the  direction  angles.  In 
other  words,  systematic  errors  on  the  long  baselines,  if  present, 
may  result  in  biased  solutions  for  the  direction  angles!  For 
data  which  is  less  than  about  6  dB  above  the  noise,  we  find  that 
both  the  RMS  deviations  and  the  dispersions  about  the  mean 
increase  as  signal  strength  decreases.  For  these  amplitudes  the 
random  component  of  the  errors  dominates  the  systematic  errors. 

At  the  station  level  of  analysis  we  have  investigated 
whether  there  are  differences  in  the  RMS  errors  for  E-W  baselines 
as  opposed  to  N-S  baselines.  We  binned  the  residuals  for  each 
individual  station  on  the  basis  of  baseline  orientation  and 
produced  polynomial  fits  to  the  data.  While  we  did  find  some 
differences  between  these  populations  at  each  station  (see  Figure 
2,  which  is  representative),  we  detected  no  repeatable  trend  from 
station  to  station.  For  each  station  we  have  plotted  the  RMS 
errors  vs.  amplitude  for  all  baselines  combined.  An  inter¬ 
comparison  of  these  plots  showed  that,  in  a  statistical  sense, 
the  error  patterns  at  all  stations  were  identical. 

The  third  level  of  analysis  which  we  performed  was  at  the 
system  level.  In  this  case  we  combined  all  the  residuals  from 
all  stations  (except  Hawkinsville)  and  binned  the  data  by  signal 
strength.  A  plot  of  the  errors  vs.  amplitude  (see  Figure  3) 
shows  that  there  are  no  significant  N-S  vs.  E-W  differences,  and 
that  the  amplitude  dependence  is  well  behaved.  In  the  Figure, 
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the  amplitude  is  expressed  in  dB  above  the  assumed  noise  floor  of 
-152  dBm.  The  curves  show  that  this  is  in  fact  a  reasonable 
estimate  of  the  overall  noise  level.  For  signals  which  are  well 
above  the  noise,  we  find  that  the  RMS  errors  are  typically  18 
degrees  and  are  nearly  independent  of  amplitude. 

The  dotted  lines  on  the  graph  show  the  best  fourth-order 
polynomial  fit  to  the  data.  The  fit  to  the  ensemble  of  all  data 
for  the  four  hour  period  investigated  is  expressed  as 

R  =  -  3.67X10'5  S4  -  1.61X10'3  S3  +  0.195  S2  -  4.44  S  +  47.2  (20) 

where  R  is  the  RMS  magnitude  of  the  residuals,  in  degrees,  and  S 
is  the  received  signal  strength  in  dB  above  the  assumed  noise 
level.  In  addition  to  examining  the  dependence  of  the  errors  on 
the  amplitude,  it  is  instructive  to  inspect  the  distribution  of 
the  phase  difference  errors  at  any  given  amplitude.  Because  the 
RMS  error  is  approximately  independent  of  signal  strength  for 
signals  greater  than  -140  dBm,  we  have  combined  these  data  to 
generate  the  histogram  in  Figure  4.  Superimposed  on  the 
histogram  is  a  Guassian  distribution  with  a  standard  deviation  of 
18  degrees.  Figure  4  demonstrates  that  the  ensemble  of  residuals 
is  nearly  normally  distributed,  although  the  skew  to  positive 
error  values  is  real  and  statistically  significant.  Similar 
curves  for  lower  signal  stregths  also  display  a  normal 
distribution  and  also  show  significant  skewing  to  positive  error 
values.  The  reason  for  this  skewing  is  not  presently  understood. 

As  part  of  our  continuing  effort  to  understand  the  nature 
and  magnitude  of  the  errors  in  the  NAVSPASUR  data,  we  are 
currently  analyzing  a  larger  data  base.  With  this  larger  data 
base  we  should  be  able  to  identify  correlations  in  the  error 
distribution  as  well  as  refine  our  error  model.  Preliminary 
results  with  the  larger  data  base  indicate  that,  for  example,  the 
errors  associated  with  the  Kickapoo  transmitter  are  larger  than 
those  associated  with  the  other  two  transmitters.  We  will 
present  the  results  of  our  further  analysis  in  a  future  report 
when  available. 

VI.  Summary  of  Results 

We  have  performed  a  nonlinear  least  squares  analysis  on  a 
four  hour  segment  of  full-Doppler  mode  NAVSPASUR  phase4  data. 
This  analysis  has  provided  direct  estimates  of  the  phase 
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difference  errors  present  in  full-Doppler  mode  data  and  the 
dependence  of  these  errors  on  received  signal  strength. 

We  find  that  for  signals  whose  amplitudes  exceed  ~  -145  dBm, 
systematic  errors,  which  are  typically  14  degrees,  dominate  the 
random  errors,  which  are  nearly  constant  over  this  amplitude 
range  and  are  typically  9  degrees.  The  systematic  errors  for  any 
given  baseline  remain  fairly  constant  over  the  time  interval  we 
have  investigated.  There  is  a  tendency  towards  larger  systematic 
errors  on  the  shortest  baselines.  However,  this  tendency  may  be 
an  artifact  of  the  least  squares  fitting  process,  since  phase 
difference  offsets  on  the  longest  baselines  result  in  biased 
solutions  to  the  direction  angles. 

At  amplitudes  below  about  -145  dBm,  the  random  component  of 
the  phase  difference  errors  dominate  the  systematic  errors,  and 
the  RMS  errors  grow  quickly  with  decreasing  amplitude.  This 
finding  is  in  agreement  with  our  assumption  that  the  effective 
noise  floor  is  approximately  -152  dBm. 

The  total  (systematic  plus  random)  phase  difference  errors 
are,  when  analyzed  as  an  ensemble,  normally  distributed  with  a 
dispersion  which  is  well  characterized  as  a  function  of  received 
signal  strength.  The  dependence  on  signal  strength  is  accurately 
modelled  by  a  fourth-order  polynomial  between  the  lowest  and 
highest  amplitudes  normally  encountered  by  the  system. 

The  random  component  of  the  phase  difference  errors  are 
comparable  in  size  to  the  uncertainties  due  to  the  6-bit 
quantization  of  the  phase  data.  However,  since  the  total  phase 
difference  errors  always  dominate  the  contribution  due  to  phase 
quantization,  the  6-bit  quantization  of  phase  data  is  not  a 
limiting  factor  in  system  performance. 


17 


REFERENCES 


Lewis,  T.  0.  ,  and  Odell,  P.  L.  1971,  Estimation  in  Linear  Models, 
(Englewood  Cliffs,  NJ :  Prentice-Hall). 

Mathews,  J.  and  Walker,  R.  L.  1970,  Mathematical  Methods  of 
Physics,  (Menlo  Park,  CA:  Benjamin  Cummings). 


18 


FIGURE  l:  Diagram  showing  the  orientation  of  a  local  station 
coordinate  system  on  the  surface  of  the  geoidal  earth.  The  h  and 
•  axes  lie  in  the  plane  of  the  NAVSPASUR  great  circle,  with  the  n 
axis  perpendicular  to  the  great  circle  plane. 
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FIGURE  2:  Plot  of  the  distribution  of  the  RMS  phase  difference 
errors  vs.  amplitude  for  the  San  Diego  receiving  station  as 
determined  from  a  least  squares  fit  to  four  hours  of  NAVSPASUR 
data.  The  dashed  lines  indicate  fourth  order  polynomial  fits  to 
the  data. 
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FIGURE  3:  Plot  of  the  RMS  phase  difference  error  distribution 
vs.  amplitude  for  all  receiving  stations  combined,  except 
Hawkinsville,  based  on  a  least  squares  fit  to  four  hours  of 
NAVSPASUR  data.  The  dashed  lines  indicate  fourth  order 
polynomial  fits  to  the  data. 
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FIGURE  4:  The  distribution  of  the  residual  phase  difference 
errors  for  data  in  the  range  -14b  to  -i^u  dBm  from  all  receiving 
stations  except  Hawkinsville .  The  data  were  derived  from  least 
squares  fits  to  four  hours  of  NAVSPASUR  data.  The  dashed  curve 
shows  a  Gaussian  distribution  of  dispersion  9.0  degrees  and  mean 
0.0  degrees. 
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PROGRAM  ALGORITHM  (PSEUDOCODE): 


ISAT  -  satellite  nuaber  . 

ISTATN  -  receiving  station  number. 

ITRAN  -  transmitter  number. 

I DOP  -  Doppler  frequency  in  Hz. 
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print  the  final  means  ar.d  RMS's.  Also,  write  the  final  residuals 
(with  headers)  to  the  individual  station  *. RSI-6  files.  If  no 
convergence,  write  the  residuals  to  the  * . DMP  file  (dump  file). 

LR  is  the  output  file  designator  for  the  residuals. 
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values  of  the  direction  angles  and  the  range. 
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INPUTS  EXPLICIT  (arguments  to  CALL  statement): 


G«t  th«  plus*  diff«c«nc»  for  PAR  -  deltaPAR. 
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